"""
    MainDriver

Main driver module for complete GSI analysis workflow orchestration,
ported from GSI's gsimain.f90 and gsisub.F90.

This module provides high-level interfaces that orchestrate the entire
data assimilation process from initialization to final output:

- Complete analysis workflow execution
- Configuration management and validation
- Module initialization and cleanup
- Error handling and diagnostics
- Integration of all GSI components

# Usage

```julia
using GSICoreAnalysis.MainDriver

# Configure analysis
config = GSIAnalysisConfig(
    analysis_type = :threedvar,
    grid_config = RegionalGridConfig(200, 150, 64),
    max_iterations = 50,
    convergence_tol = 1e-6
)

# Run complete analysis
result = run_gsi_analysis(config, 
                         background_file="wrfout_d01.nc",
                         observation_file="prepbufr.dat", 
                         output_file="analysis.nc")
```
"""
module MainDriver

using Printf
using Dates

# Import all GSI components
using ..GSICoreAnalysis
using ..GridOperations
using ..DataIO  
using ..BackgroundError
using ..ObservationOperators
using ..CostFunctions
using ..Minimization

# Export main functions
export GSIAnalysisConfig, GSIAnalysisResult
export run_gsi_analysis, validate_configuration
export initialize_gsi_system, finalize_gsi_system

"""
    GSIAnalysisConfig

Complete configuration for GSI analysis run.

# Fields
- `analysis_type::Symbol`: Type of analysis (:threedvar, :fourdvar, :hybrid_threedvar, :hybrid_fourdvar)
- `grid_config::AbstractGridConfig`: Grid configuration
- `background_files::Vector{String}`: Background field files
- `observation_files::Vector{String}`: Observation data files  
- `output_file::String`: Analysis output file
- `diagnostic_file::String`: Diagnostic output file
- `solver_config::Dict{Symbol,Any}`: Minimization solver configuration
- `quality_control::Bool`: Enable observation quality control
- `hybrid_ensemble::Bool`: Enable hybrid ensemble-variational method
- `ensemble_size::Int`: Number of ensemble members (if hybrid)
- `localization_radius::Float64`: Localization radius (km)
- `inflation_factor::Float64`: Ensemble inflation factor
"""
struct GSIAnalysisConfig
    analysis_type::Symbol
    grid_config::AbstractGridConfig
    background_files::Vector{String}
    observation_files::Vector{String}
    output_file::String
    diagnostic_file::String
    solver_config::Dict{Symbol,Any}
    quality_control::Bool
    hybrid_ensemble::Bool
    ensemble_size::Int
    localization_radius::Float64
    inflation_factor::Float64
end

# Outer constructor with defaults
function GSIAnalysisConfig(
    grid_config::AbstractGridConfig,
    background_files::Vector{String},
    observation_files::Vector{String};
    analysis_type::Symbol = :threedvar,
    output_file::String = "analysis.nc",
    diagnostic_file::String = "gsi_diagnostics.nc",
    solver_config::Dict{Symbol,Any} = Dict(:max_iterations => 50, :tolerance => 1e-6),
    quality_control::Bool = true,
    hybrid_ensemble::Bool = false,
    ensemble_size::Int = 0,
    localization_radius::Float64 = 500.0,
    inflation_factor::Float64 = 1.05)
    
    # Validation
    @assert analysis_type in [:threedvar, :fourdvar, :hybrid_threedvar, :hybrid_fourdvar] "Invalid analysis type"
    @assert !isempty(background_files) "Must specify background files"
    @assert !isempty(observation_files) "Must specify observation files"
    
    if hybrid_ensemble
        @assert ensemble_size > 0 "Ensemble size must be positive for hybrid methods"
    end
    
    return GSIAnalysisConfig(
        analysis_type, grid_config, background_files, observation_files,
        output_file, diagnostic_file, solver_config, quality_control,
        hybrid_ensemble, ensemble_size, localization_radius, inflation_factor)
end

"""
    GSIAnalysisResult

Result structure for complete GSI analysis.

# Fields
- `success::Bool`: Whether analysis completed successfully
- `final_cost::Float64`: Final cost function value
- `iterations::Int`: Total iterations performed
- `analysis_time::DateTime`: Time when analysis was performed
- `convergence_info::Dict`: Detailed convergence information
- `diagnostic_info::Dict`: Diagnostic information
- `timing_info::Dict`: Performance timing information
"""
struct GSIAnalysisResult
    success::Bool
    final_cost::Float64
    iterations::Int
    analysis_time::DateTime
    convergence_info::Dict{String,Any}
    diagnostic_info::Dict{String,Any}
    timing_info::Dict{String,Any}
end

"""
    run_gsi_analysis(config::GSIAnalysisConfig) -> GSIAnalysisResult

Execute complete GSI analysis workflow.

This is the main entry point that orchestrates the entire data assimilation process:
1. Initialize GSI system components
2. Read background fields and observations  
3. Setup grid operations and domain decomposition
4. Initialize background error covariance
5. Setup observation operators
6. Construct and minimize cost function
7. Write analysis output and diagnostics
8. Cleanup and finalize

# Arguments
- `config`: Complete GSI analysis configuration

# Returns
- Analysis result with convergence and diagnostic information

# Example
```
config = GSIAnalysisConfig(
    grid_config = RegionalGridConfig(200, 150, 64),
    background_files = ["background.nc"],
    observation_files = ["obs.dat"]
)
analysis_result = run_gsi_analysis(config)
println("Analysis ", analysis_result.success ? "succeeded" : "failed")
```
"""
function run_gsi_analysis(config::GSIAnalysisConfig)
    analysis_start_time = now()
    timing_info = Dict{String,Any}()
    
    println("="^60)
    println("GSI Analysis Driver - Julia Implementation")
    println("Analysis type: $(config.analysis_type)")
    println("Start time: $(analysis_start_time)")
    println("="^60)
    
    try
        # Phase 1: System Initialization
        println("\n--- Phase 1: System Initialization ---")
        t1 = time()
        
        validate_configuration(config)
        gsi_components = initialize_gsi_system(config)
        
        timing_info["initialization"] = time() - t1
        println("Initialization completed in $(timing_info["initialization"]:.2f) seconds")
        
        # Phase 2: Data Input
        println("\n--- Phase 2: Data Input ---") 
        t2 = time()
        
        # Read background fields
        background_data = read_background_data(config, gsi_components[:grid])
        
        # Read observations
        observation_data = read_observation_data(config)
        
        # Apply quality control
        if config.quality_control
            observation_data = apply_observation_qc(observation_data)
        end
        
        timing_info["data_input"] = time() - t2
        println("Data input completed in $(timing_info["data_input"]:.2f) seconds")
        println("Background fields loaded: $(length(config.background_files))")
        println("Observations loaded: $(length(observation_data.obs_values))")
        
        # Phase 3: Setup Analysis System
        println("\n--- Phase 3: Analysis System Setup ---")
        t3 = time()
        
        # Setup background error covariance
        background_error = setup_background_error(config, gsi_components[:grid])
        
        # Setup observation operators
        obs_operators = setup_observation_operators(config, gsi_components[:grid], observation_data)
        
        # Create cost function
        cost_function = create_cost_function(config, background_data, observation_data, 
                                           background_error, obs_operators)
        
        timing_info["setup"] = time() - t3
        println("Analysis setup completed in $(timing_info["setup"]:.2f) seconds")
        
        # Phase 4: Analysis (Minimization)
        println("\n--- Phase 4: Analysis Minimization ---")
        t4 = time()
        
        # Create initial guess (control vector)
        initial_guess = create_initial_guess(config, background_data)
        
        # Setup solver
        solver = create_analysis_solver(config)
        
        # Run minimization
        println("Starting $(typeof(solver)) minimization...")
        min_result = minimize_cost_function(cost_function, initial_guess, solver)
        
        timing_info["analysis"] = time() - t4
        println("Analysis completed in $(timing_info["analysis"]:.2f) seconds")
        println("Converged: $(min_result.converged)")
        println("Iterations: $(min_result.iterations)")
        println("Final cost: $(min_result.final_cost)")
        
        # Phase 5: Output and Diagnostics
        println("\n--- Phase 5: Output Generation ---")
        t5 = time()
        
        # Convert solution to analysis fields
        analysis_fields = solution_to_analysis_fields(min_result.solution, background_data, 
                                                     gsi_components[:grid])
        
        # Write analysis output
        write_analysis_output(config, analysis_fields)
        
        # Write diagnostics
        diagnostic_info = create_diagnostic_info(min_result, observation_data, cost_function)
        write_diagnostic_output(config, diagnostic_info)
        
        timing_info["output"] = time() - t5
        println("Output generation completed in $(timing_info["output"]:.2f) seconds")
        
        # Phase 6: Cleanup
        println("\n--- Phase 6: Cleanup ---")
        finalize_gsi_system(gsi_components)
        
        # Create result
        total_time = time() - t1
        timing_info["total"] = total_time
        
        convergence_info = Dict{String,Any}(
            "converged" => min_result.converged,
            "iterations" => min_result.iterations,
            "final_cost" => min_result.final_cost,
            "final_gradient_norm" => min_result.final_gradient_norm,
            "cost_history" => min_result.cost_history
        )
        
        result = GSIAnalysisResult(
            true,  # success
            min_result.final_cost,
            min_result.iterations,
            analysis_start_time,
            convergence_info,
            diagnostic_info,
            timing_info
        )
        
        println("\n" * "="^60)
        println("GSI Analysis COMPLETED SUCCESSFULLY")
        println("Total runtime: $(total_time:.2f) seconds")
        println("Final cost: $(min_result.final_cost)")
        println("Analysis file: $(config.output_file)")
        println("Diagnostic file: $(config.diagnostic_file)")
        println("="^60)
        
        return result
        
    catch e
        println("\n" * "!"^60)
        println("GSI Analysis FAILED")
        println("Error: $e")
        println("!"^60)
        
        # Return failure result
        return GSIAnalysisResult(
            false,  # success
            Inf,    # final_cost
            0,      # iterations
            analysis_start_time,
            Dict{String,Any}("error" => string(e)),
            Dict{String,Any}(),
            timing_info
        )
    end
end

# Implementation functions

"""
    validate_configuration(config::GSIAnalysisConfig) -> Nothing

Validate analysis configuration and check for required files.
"""
function validate_configuration(config::GSIAnalysisConfig)
    println("Validating configuration...")
    
    # Check background files exist
    for file in config.background_files
        if !isfile(file)
            error("Background file not found: $file")
        end
    end
    
    # Check observation files exist
    for file in config.observation_files
        if !isfile(file)
            error("Observation file not found: $file")
        end
    end
    
    # Validate solver configuration
    if !haskey(config.solver_config, :max_iterations)
        error("Solver configuration missing max_iterations")
    end
    
    if !haskey(config.solver_config, :tolerance)
        error("Solver configuration missing tolerance")  
    end
    
    println("Configuration validation passed")
end

"""
    initialize_gsi_system(config::GSIAnalysisConfig) -> Dict{Symbol,Any}

Initialize all GSI system components.
"""
function initialize_gsi_system(config::GSIAnalysisConfig)
    components = Dict{Symbol,Any}()
    
    println("Initializing grid operations...")
    grid = initialize_grid(config.grid_config)
    components[:grid] = grid
    
    # Initialize domain decomposition if needed
    if isa(config.grid_config, RegionalGridConfig)
        println("Setting up domain decomposition...")
        # For now, single processor
        components[:decomposition] = nothing
    end
    
    println("GSI system initialization completed")
    return components
end

"""
    finalize_gsi_system(components::Dict{Symbol,Any}) -> Nothing

Cleanup and finalize GSI system components.
"""
function finalize_gsi_system(components::Dict{Symbol,Any})
    println("Finalizing GSI system...")
    # Cleanup operations would go here
    println("GSI system finalized")
end

# Placeholder implementations for component functions

function read_background_data(config::GSIAnalysisConfig, grid::GridOperations.AbstractGrid)
    println("Reading background data from $(length(config.background_files)) files...")
    # Use DataIO module
    return read_background_fields(config.background_files[1])
end

function read_observation_data(config::GSIAnalysisConfig)
    println("Reading observation data from $(length(config.observation_files)) files...")
    # Use DataIO module  
    return read_observations(config.observation_files[1])
end

function apply_observation_qc(obs_data)
    println("Applying observation quality control...")
    # Use ObservationOperators module
    return obs_data  # Simplified
end

function setup_background_error(config::GSIAnalysisConfig, grid::GridOperations.AbstractGrid)
    println("Setting up background error covariance...")
    # Use BackgroundError module
    return "background_error_placeholder"
end

function setup_observation_operators(config::GSIAnalysisConfig, grid::GridOperations.AbstractGrid, obs_data)
    println("Setting up observation operators...")
    # Use ObservationOperators module
    return [SurfacePressureOperator(grid)]
end

function create_cost_function(config, bg_data, obs_data, bg_error, obs_ops)
    println("Creating cost function...")
    # Use CostFunctions module
    return "cost_function_placeholder"
end

function create_initial_guess(config::GSIAnalysisConfig, bg_data)
    println("Creating initial guess (control vector)...")
    # Use ControlVectors module
    return "initial_guess_placeholder"
end

function create_analysis_solver(config::GSIAnalysisConfig)
    println("Creating analysis solver: $(get(config.solver_config, :type, :pcg))")
    
    # Create appropriate solver based on analysis type
    if config.analysis_type in [:threedvar, :hybrid_threedvar]
        if config.hybrid_ensemble
            return GlobalAnalysisSolver(config.grid_config; 
                                      max_outer_iterations=1,
                                      hybrid_ensemble=true)
        else
            return PCGSolver(config.grid_config; 
                           max_iterations=config.solver_config[:max_iterations],
                           tolerance=config.solver_config[:tolerance])
        end
    elseif config.analysis_type in [:fourdvar, :hybrid_fourdvar]
        return GlobalAnalysisSolver(config.grid_config;
                                  max_outer_iterations=3,
                                  four_dvar=true,
                                  hybrid_ensemble=config.hybrid_ensemble)
    else
        error("Unknown analysis type: $(config.analysis_type)")
    end
end

function solution_to_analysis_fields(solution, bg_data, grid::GridOperations.AbstractGrid)
    println("Converting solution to analysis fields...")
    # Convert control vector solution to physical fields
    return bg_data  # Simplified - would apply increments
end

function write_analysis_output(config::GSIAnalysisConfig, analysis_fields)
    println("Writing analysis output to $(config.output_file)...")
    # Use DataIO module
    # write_analysis_fields(config.output_file, analysis_fields)
end

function create_diagnostic_info(min_result, obs_data, cost_function)
    return Dict{String,Any}(
        "observations_used" => length(obs_data.obs_values),
        "observations_passed_qc" => length(obs_data.obs_values),
        "cost_reduction" => length(min_result.cost_history) > 1 ? 
                           min_result.cost_history[1] - min_result.cost_history[end] : 0.0
    )
end

function write_diagnostic_output(config::GSIAnalysisConfig, diagnostic_info)
    println("Writing diagnostic output to $(config.diagnostic_file)...")
    # Create and write diagnostic file
end

end # module MainDriver